/**
  SUNNY Project, Anyang Normal University, IMP-CAS
  \class TAMBSDManager
  \brief A class to generate many-body basis and manage TAMBSDList objects.
  \author SUN Yazhou
  \since 2020/02/01
  \date Last modified: 2022/11/15 by SUN Yazhou
	\copyright SUNNY project, Anyang Normal University, IMP-CAS
*/

#include <iostream>
#include <algorithm>
#include <cstring>
#include <vector>
#include <unistd.h>
#include "config.h" // contains user switch of definition of ROOT and OpenCV
#include "TAMBSDManager.h"
#include "TAMBSDList.h"
#include "TAMBSD.h"
#include "TASPStateManager.h"
#include "TAException.h"
#include "TAYaml.h"

using std::cout;
using std::endl;
using std::vector;

TAMBSDManager *TAMBSDManager::kInstance = nullptr;

TAMBSDManager::TAMBSDManager() : fMBSDListM(nullptr),
    fNParticle(0), fNSPState(0), f2M(-9999){
  fCfgFile = "";
  fBothParity = false; fIsPair = false;
  fParity = true; // default is true -> +
} // end the constructor

TAMBSDManager *TAMBSDManager::Instance(){
  if(!kInstance) kInstance = new TAMBSDManager();
  return kInstance;
}

TAMBSDManager::~TAMBSDManager(){}

bool TAMBSDManager::MBSDFilter(ull bit) const{
  TAMBSD m(bit);
  // return true; // no restrictions -- DEBUG
  // M-sheme and T-scheme: M and T selection
  const bool mt = (m.Get2T() == f2T) && (m.Get2M() == f2M); // test total M and T
  // parity selection
  const bool p = fBothParity ? true : (fParity == m.Parity());
  // if Pairing is required, m.IsPaired() must be true
  // or the value of m.IsPaired() is irrelevant
  bool pair = fIsPair ? m.IsPaired() : true;
  return
    pair &&
    mt &&
    p &&
    true;
} // end member function MBSDFilter

// skip spaces and tabs, return subscript of the valid char
inline int skipCrap(const char *s){
	int tmp = 0;
	while(1){
		char c = s[tmp++];
		if(' ' != c && '\t' != c) break;
	}
	return tmp - 1;
}
// assign nparticle, nspstate, 2*M and the spstate space
void TAMBSDManager::LoadConfigFile(const string &file){
  if(access(file.c_str(), F_OK)) TAException::Error("TAMBSDManager",
     "Read: file %s non-existent", file.c_str());
  fCfgFile = file;
  TAYaml pa(file);
  fNParticle = pa["N"].as<int>();
  f2M = pa["2M"].as<int>();
  f2T = pa["2T"].as<int>();
  if(pa["pi"]) fParity = pa["pi"].as<bool>();
  else fBothParity = true;
  if(pa["IsPair"]) fIsPair = pa["IsPair"].as<bool>();
} // end member function LoadConfigFile


void TAMBSDManager::GenerateManyBodySD(){
  if(fMBSDListM->GetNBasis()) return; // called already

  // generate SP state, MB state, and M-scheme MB state list //
  TASPStateManager *spStateManager = TASPStateManager::Instance();
  spStateManager->LoadSPListFile(fCfgFile);
  fNSPState = spStateManager->GetNSPState();
  if(fNParticle > fNSPState) TAException::Error("TAMBSDManager",
    "GenerateManyBodySD: nParticle > nSPState.");

  /////////// odometer method to generate many-body basis /////////////
  if(0 == fNParticle) TAException::Error("TAMBSDManager",
    "GenerateManyBodySD: Number of particles not assigned.");
  vector<int> spStateVec; spStateVec.resize(fNParticle);
  // the first MBSD configuration
  for(int i = 0; i < fNParticle; i++) spStateVec[i] = i;
  fMBSDListM->Add(spStateVec.data(), fNParticle); ///< Fill mbsd in a TTree obj

  // generate the many-body slater determinants //
  // fNParticle sp-states are stored in each int of array spStateVec. The sp-states are
  // labeled from 0, and all the states are distinct, e.g. [0,1,2,3], [0,1,2,4], ...
  // the odometer method starts from the last sp-state, increments it to the maximum
  // (maxSP), and then recover the following sp-states. e.g. [3,4,5,6]->[4,5,6,7]
  // for 4 particles in 8 sp-states.
  while(fNSPState - fNParticle != spStateVec[0]){
    int i = fNParticle - 1, maxSP = fNSPState - 1; // max state of the last particle
    // find the bit that is not full (overflow) in array spStateVec backwards
    while(i >= 0 && spStateVec[i--] == maxSP--);
    spStateVec[++i]++; // increment the found unfull bit
    // recover the following (rightward) sp-states to reinitiate incrementation //
    while(i < fNParticle - 1){ spStateVec[i + 1] = spStateVec[i] + 1; i++; }
    fMBSDListM->Add(spStateVec.data(), fNParticle); ///< register the mbsd
  } // end while
  ////////////////// END of the odometer algorithm /////////////////////
} // end member function GenerateManyBodySD

// generate the M-scheme many-body state basis
void TAMBSDManager::MSchemeGo(){
  if(fMBSDListM) return; // already called

  if(-9999 == f2M) TAException::Error("TAMBSDManager",
    "GenerateManyBodySD: f2M for M-scheme not assigned.");
  fMBSDListM = new TAMBSDList();
  GenerateManyBodySD(); // Generate all the many-body basis, plus filter

  if(0 == fMBSDListM->GetNBasis()) TAException::Warn("TAMBSDManager",
    "MSchemeGo: fMBSDListM is empty in the end.");
#ifdef _ROOT_
  fMBSDListM->Save(); // save the mbsd TTree obj
#endif
} // end member function MSchemeGo

TAMBSDList *TAMBSDManager::GetMBSDListM(){
  if(!fMBSDListM) TAException::Warn("TAMBSDManager",
    "GetMBSDListM: fMBSDListM nullptr. MSchemeGo() not called yet?");
  return fMBSDListM;
} // end member function GetMBSDListM

// display the geneated many-body SD
// opt: 1->Print; 2->PrintInBit; 3->Print both; otherwise->brief;
void TAMBSDManager::Print(int opt){
  cout << "\033[33;1m--- MBSD Manager Info ---\033[0m" << endl;
  cout << "Number of particles: " << fNParticle << endl;
  if(fBothParity) cout << "Both + and - parities are allowed." << endl;
  else cout << "Parity: " << ((fParity ? '+' : '-')) << endl;
  cout << "2M for M-scheme: " << f2M << endl;
  cout << "2T for T-scheme: " << f2T << endl;

	const int n = fMBSDListM->GetNBasis();
	if(n <= 0) return;
	cout << "\033[33;1mTotally " << n << " many-body Slater determinants ";
  cout << "accepted in the list.\033[0m" << endl;

  if(1 == opt) fMBSDListM->Print();
  if(2 == opt) fMBSDListM->PrintInBit();
  if(3 == opt){
    fMBSDListM->Print();
    fMBSDListM->PrintInBit();
  } // end if
} // end member function Print
